%% 对带有边界的一元数值函数求其泰勒展开式
%% 仿照：taylor(f,x,a,'order',n)
%[A,coefficient]=Taylor_series(f,x0,n,xmin,dx,xmax)
function [A,coefficient]=Taylor_series(f,x0,n,xmin,dx,xmax)
%% 函数原始调试
% n=6;%泰勒级数阶数
% x0=0;%求级数的点位置
% xmin=-10;xmax=10;dx=1e-4;
% x=xmin:dx:xmax;
% f=sin(x);
%% 函数说明：
% 输入量：f为函数信号，为数组；x0为泰勒展开的自变量取值点；n为泰勒级数阶数；
% xmin为自变量取值最小值；xmax为自变量取值最大值，dx为取值间隔。
% 输出量：A为各阶泰勒级数，A第i行代表第i阶泰勒展开式；coefficient为各阶泰勒
% 级数系数，第i项代表第i阶泰勒展开式系数
%% 函数注意：
% n值不能过大，5左右效果较优；
% 输入f量不能过大。
%% 构造函数
coefficient=zeros(1,n);% 泰勒级数各阶系数
coefficient(1)=f(round((x0-xmin)/dx)+1);
a0=diff(f)/dx;% 利用差分求近似导数
k1=1;k2=1;% 对于不同次导数，处理方法不同
for i=2:n
    a1=length(a0);
    k0=mod(i,2);
    if(k0==0)% 当导数次数为奇数次时，求该点处两端差分均值
        if(round((x0-xmin)/dx)-k1+1>=0&&round((x0-xmin)/dx)-k1+2<=a1)
            coefficient(i)=(a0(round((x0-xmin)/dx)-k1+1)+...
                a0(round((x0-xmin)/dx)-k1+2))/2;
        elseif(round((x0-xmin)/dx)-k1+1<0)
                coefficient(i)=a0(1);
        else
                coefficient(i)=a0(end);
        end
        k1=k1+1;
    else% 当导数次数为偶数次时，求该点差分
        if(round((x0-xmin)/dx)-k2>=0&&round((x0-xmin)/dx)-k2<=a1)
            coefficient(i)=a0(round((x0-xmin)/dx)-k2);
        elseif(round((x0-xmin)/dx)-k2<0)
            coefficient(i)=a0(1);
        else
            coefficient(i)=a0(end);
        end
        k2=k2+1;
    end
    a0=diff(a0)/dx;
end
for i=1:n
    coefficient(i)=coefficient(i)/factorial(i-1);
end
xm=xmin-abs(x0):dx:xmax+abs(x0);
length_xm=length(xm);
A=zeros(n,length_xm);
for i=1:n
    A(i,:)=coefficient(i)*(xm-x0).^(i-1);
end
end